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Abstract 

Background: Understanding the normal temporal variation in the human microbiome is critical to developing 
treatments for putative microbiome-related afflictions such as obesity, Crohn's disease, inflammatory bowel disease 
and malnutrition. Sequencing and computational technologies, however, have been a limiting factor in performing 
dense time series analysis of the human microbiome. Here, we present the largest human microbiota time series 
analysis to date, covering two individuals at four body sites over 396 timepoints. 

Results: We find that despite stable differences between body sites and individuals, there is pronounced variability 
in an individual's microbiota across months, weeks and even days. Additionally, only a small fraction of the total 
taxa found within a single body site appear to be present across all time points, suggesting that no core temporal 
microbiome exists at high abundance (although some microbes may be present but drop below the detection 
threshold). Many more taxa appear to be persistent but non-permanent community members. 

Conclusions: DNA sequencing and computational advances described here provide the ability to go beyond 
infrequent snapshots of our human-associated microbial ecology to high-resolution assessments of temporal 
variations over protracted periods, within and between body habitats and individuals. This capacity will allow us to 
define normal variation and pathologic states, and assess responses to therapeutic interventions. 



Background 

As more attention is paid to viewing ourselves as a 
supraorganism, comprising interacting microbial and 
human cellular and genetic components, it is apparent 
that much more precise understanding is needed of 
what constitutes normal temporal variations in our 
microbial community structures and functions. Variation 
in the human microbiome within and between our var- 
ious body habitats, lifecycle stages, and cultural settings 
is largely unexplored. High-resolution time series studies 
provide a foundation for discriminating between 'nor- 
mal' perturbations and pathologic states, and between 
organisms that are simply passing through a body habi- 
tat or are entrenched residents of an ecosystem. Simi- 
larly, these types of studies are needed to understand 
the immigration and emigration patterns of microbes 
between our body sites, between cohabitating indivi- 
duals, and between ourselves and the myriad of environ- 
ments we contact on a daily basis [1-3]. 
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The densest human microbiome time series reported 
to date studied the response of distal gut microbial com- 
munities to the antibiotic ciprofloxacin across three 
individuals, with sampling intervals varied from daily to 
weekly. Eighteen timepoints were collected per indivi- 
dual in an initial study [4], and between 52 and 56 time- 
points were collected per subject in a follow-up study 
[5]. In both studies, rapid decreases in alpha diversity 
and a characteristic shift in community composition 
were observed in association with antibiotic therapy, fol- 
lowed by a rapid post-antibiotic increase in diversity as 
the gut community returned to a state similar (but not 
identical) to the pre-treatment state. Another dense 
human microbiome time series studied 60 fecal samples 
over the first 2.5 years of life for a single infant, and illu- 
strated the successional pattern in the human gut 
microbiota as a rapidly changing community that devel- 
ops over time into a community characteristic of that 
found in the adult human gut [2]. 

The availability of radically cheaper sequencing and 
analysis, described here, eases constraints on the num- 
ber of timepoints and body sites that can be compared 
in a single study, and paves the way to a finer-grained 
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understanding of how the human microbiota changes in 
different body habitats (including much better estimates 
of the relative variability of different body sites within a 
subject). Understanding this intrinsic variability will be 
crucial for performing power calculations to test 
whether antibiotics, probiotics, or other drugs actually 
affect the microbiome. 

In this study, two healthy subjects, one male (M3) and 
one female (F4), one of whom (M3) participated in an 
earlier survey [6], were sampled daily at three body sites 
(gut (feces), mouth, and skin (left and right palms)), for 
15 months (M3) and for 6 months (F4) using an institu- 
tional review board-approved protocol Variable region 
4 (V4) of 16S rRNA genes present in each community 
sample were amplified by PCR and subjected to multi- 
plex sequencing on an Illumina Genome Analyzer IIx 
(GA-IIx; average read length after quality trimming, 123 
± 17 (standard deviation (SD)) nucleotides; 32,266 ± 
19,723 (SD) reads per sample; n = 1,422 (M3 samples); 
n - 531 (F4 samples); average interval between sam- 
pling, 1.12 days). To control for differences across 
sequencing platforms and primer pairs, 331 of these 
samples had variable region 2 (V2) sequenced on 454 
(average read length after quality filtering, 228 ±11 
(SD) nucleotides; 1,072 ± 375 (SD) reads per sample; 
n - 171 (M3 samples); n - 160 (F4 samples)). This 
study thus provides a key counterpoint to recent studies 
of one or a few subjects at tens of timepoints [2-5,7], or 
studies that examine hundreds of individuals but only at 
one or a few timepoints [6,8-12]. 

Results and discussion 

Stable differences in microbial communities between 
body sites over time 

When the samples from Costello et al. [6] and the cur- 
rent study are compared directly, the samples cluster by 
body habitat, showing excellent concordance between 
the studies (Figure la in [6]; Figure lb in the current 
study) despite differences in sequencing technology (454 
and Illumina GA-IIx, respectively), mean read length 
(229 ± 16 (SD) nucleotides and 123 ± 17 (SD) nucleo- 
tides, respectively), and region of the 16S sequenced (V2 
and V4, respectively). The UniFrac distances between 
the 331 time series samples that were sequenced on 
both Illumina and 454 were significantly correlated, as 
determined by Procrustes analysis of unweighted Uni- 
Frac principal coordinate matrices (M2 = 0.161; Monte 
Carlo P < 0.001; Additional file 1) and Pearson correla- 
tion of UniFrac distances for pairs of samples (r = 0.91; 
P < 0.001). As observed by Costello et al. [6], gut, oral, 
and skin bacterial communities were found to be com- 
positionally distinct based on principal coordinates ana- 
lysis of unweighted UniFrac distances between 
communities (UniFrac measures community similarity 



Figure 1 Principal coordinates analysis of unweighted UniFrac 
distances between samples, (a) Costello et al. [6] samples, (b) M3, 
F4 time series samples, (c) M3, F4 time series, PC1 versus time 
(days). Panels (a,b) and (c) show two independent principal 
coordinates analyses. To compare the Costello et al. 454 data (a) 
with the time series Illumina data (b), these data were generated in 
a single principal coordinates analysis of UniFrac distances at 500 
sequences per sample. Panel (c) does not contain the 454 data, so 
makes use of the increased sampling depth possible on Illumina 
(evenly sampled to 5,000 sequences per sample for UniFrac 
calculations). 

L J 

based on the degree to which they share branch length 
on a phylogenetic tree). The long-term time series 
shows for the first time that this body-site differentiation 
is highly stable over greater than one year (Figure lc), 
but dynamic within sites over time (Additional files 2, 3, 
4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15). 

Minimal evidence for a temporal core microbiome 
between or within body sites 

While overall compositional differences between body 
sites and individuals were relatively stable, our data also 
suggest a surprisingly small temporal 'core human 
microbiota' within an individual's body sites (Figure 2) 
when we define the 'core' as those species-level phylo- 
types in a given body habitat that were observed across 
all sampling events. These data suggest a minimal core 
microbiome across time, where the size of the core 
decreases as: mouth > gut > right palm « left palm > 
across body sites within an individual > across body 
sites and individuals. 

At this depth of sequencing, many more OTUs are 
either persistent community members, which appear in 
a given body habitat and remain for an extended period 
of time but are not present consistently enough to be 
considered core members, or transient community 
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Minimum number of samples to be considered part of core 

Figure 2 Temporal core microbiome. Fraction of species-level 
operational taxonomic units (OTUs) composing the core microbiota 
by number of samples in which an OTU must be present to be 
considered part of the core. 



members, which appear in a body habitat and disappear 
soon after (Figure 3). The taxa composing these persis- 
tent and transient categories are significantly different 
(for example, M3 gut: G indep , 84.78; P = 1.80 x 10" 14 ). In 
M3 gut, both the persistent and transient communities 
are dominated by Clostridia, Bacteroidia, and to a lesser 
extent Erysipelotrichi. The persistent community is, 
however, also composed of Betaproteobacteria and Del- 
taproteobacteria, while the transient community is com- 
posed of Actinobacteria, Gammaproteobacteria, 
Epsilonproteobacteria, and Verrucomicrobiae. Taxo- 
nomic summaries of the persistent and transient groups 
for all body sites from both individuals are presented in 
Additional file 16. 

We applied several techniques to control for the pos- 
sibility of persistent groups being mistaken for transient 
groups if they occasionally fell below the detection 
threshold. First, in computing the maximum number of 
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Figure 3 Community membership. Community membership summary for all OTUs in (a) M3 gut, (b) F4 gut, (c) M3 tongue, (d) F4 tongue, (e) 
M3 left palm, (f) F4 left palm, (g) M3 right palm, and (h) F4 right palm. Points are OTUs colored by their median relative abundance computed 
over all samples where they occur, and pie charts summarize the class-level taxa observed as persistent and transient OTUs. 
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consecutive observations for an OTU, we counted a sin- 
gle zero-count for an OTU as not interrupting a run of 
consecutive observations, provided that both adjacent 
timepoints achieved non-zero counts for that OTU. Sec- 
ond, we re-sequenced 331 of the time series samples on 
the 454 platform, and recomputed the persistent and 
transient taxa summaries. We found that, for both indi- 
viduals, the composition of the persistent gut and oral 
communities were not significantly different between 
454 and Illumina despite a nearly 40-fold difference in 
sequencing depth. However, the persistent palm com- 
munities were significantly different in both cases. A 
possible confounding factor in this comparison is primer 
bias, as the V2 region was sequenced on 454. We there- 
fore performed 1,000 jackknife iterations of this analysis 
by subsampling the Illumina data to 5,000 sequences 
per sample and recomputing the composition of the 
persistent group for each individual and body habitat. In 
this analysis, seven of the eight individual/body site 
pairs were never significantly different from the persis- 
tent community composition on the full Illumina data 
set. The one exception was the F4 gut community, 
which was significantly different in approximately half of 
the iterations. When re-sampled at a depth of 10,000 
sequences per sample, the persistent F4 gut community 
never achieved significant difference from that deter- 
mined on the full Illumina data set. The full results of 
these analyses are presented in Additional file 17. The 
designation of groups as 'persistent' was thus highly 
reproducible across sequencing platforms and ampli- 
cons, although it is still possible that sequencing error 
or very low abundance taxa occasionally falling below 
the detection threshold could result in underestimated 
size of the persistent group. 

Temporally dynamic microbial communities and 
correlations between body sites 

Differences in UniFrac distances in the left and right 
palms in adjacent timepoints of both individuals were 
significantly correlated (Pearson correlation for M3, r = 
0.69, P = 2.07 x 10~ 46 ; for F4, r = 0.64, P = 3.77 x 10" 
16 ), possibly due to equilibration of microbial commu- 
nities across palms by physical contact. We did not see 
correlations between other body sites. While the magni- 
tude and direction of change in phylogenetic dissimilar- 
ity between adjacent timepoints were correlated between 
the palm sites, the species-level microbial taxa present 
on each hand were not significantly correlated, confirm- 
ing previous observations that, at a single timepoint, the 
left and right hand of a single individual may share rela- 
tively few OTUs [13]. 

For single body sites, within-subject distances were 
lower than between-subject distances, suggesting a 
stable pattern, consistent across body sites, in between- 



subject dissimilarities across time. For example, the 
between-subject fecal sample distances were significantly 
higher than the within M3 fecal sample distances (t = 
15.52; P < 0.001; one-tailed, two sample £-test) and the 
within F4 fecal sample distances (t = 33.45; P < 0.001; 
one- tailed, two sample £-test). 

The microbial community dynamics are especially 
apparent in principal coordinates analysis animations 
(Additional files 2, 3, 4, 5, 6, 7) where the sample types 
(subject, body site combinations) are represented as 
moving traces against a background of the Costello et 
al. data. The traces give an immediate picture of the 
variability within each body site, the relative distinctive- 
ness of the sites and the subjects, and the relative speed 
of change in each site. Blooms of certain taxa contribute 
to the within site dissimilarities over time, as with the 
waning and waxing of the relative Proteobacteria abun- 
dance in the gut of both M3 and F4 (Additional files 8 
and 9, Phylum panel (teal)). Similar patterns are appar- 
ent across all taxonomic levels and body sites (Addi- 
tional files 8, 9, 10, 11, 12, 13, 14, 15). Such 
visualizations of microbial community dynamics will add 
another dimension to long-term studies of variable clini- 
cal states, such as inflammatory bowel disease or drug 
treatments, and changes in diet or lifestyle. 

Conclusions 

Dense, deeply sequenced microbiome time series studies 
have been limited by the cost of sequencing, and the 
computational power necessary to analyze such studies 
has been expensive. Several technical advances made 
this study possible. First, using the Illumina sequencer 
and the protocol described in [14] reduced sequencing 
cost per sample ($11 USD) to below the cost associated 
with DNA extraction and PCR ($13). Second, cloud 
computing using Amazon Web Services (AWS) allowed 
us to perform the bioinformatics analysis for $200 using 
open-source and freely available software running on 
commodity services. This included the creation of a vir- 
tual cluster on AWS composed of 20 eight-core systems 
with 68 GB of RAM each (that is, 20 Amazon 
m2.4xlarge instances) to cluster approximately 69 mil- 
lion sequences into OTUs. This step alone would have 
required 20 days of computation on a desktop machine 
with 8 GB of RAM, but was performed on this virtual 
cluster in 3 hours for $120. 

Figure 2 illustrates the sensitivity of conclusions about 
the core microbiome to the definition of the core micro- 
biome: if an OTU is considered here to be a core mem- 
ber only if it is present in all 130 samples, as compared 
to 120 of the samples, the size of the core drops to 
approximately 5% of OTUs from about 10% of OTUs. 
This rapid change in slope beyond 120 samples on the 
x-axis does not occur when considering all M3 gut or 
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oral samples rather than only 130 (data not shown), sug- 
gesting that this result likely stems from a few non- 
representative samples, rather than representing an 
interesting biological result. The definition of the core 
microbiome is thus a crucial consideration, because 
depth of sampling, PCR error, sequencing error, analysis 
techniques (such as OTU picking method and similarity 
threshold), among other factors, will affect whether an 
OTU will be observed in a given sample. OTU picking 
similarity threshold, in particular, will have a large affect 
on the size of the observed core microbiome. As this 
threshold is decreased, OTUs represent larger taxo- 
nomic groups, and the observed core will increase in 
size correspondingly. For example, there is clearly a 
temporal core at the phylum level in the fecal samples 
studied here composed of Bacteriodetes and Firmicutes 
(Additional files 8 and 9). For this reason our results do 
not necessarily contradict previous Denaturing Gradient 
Gel Electrophoresis (DGGE)-based work in this area 
[4,15-17], which is suggestive of a core gut microbiome 
but at a coarser phylogenetic resolution and more lim- 
ited dynamic range. Sampling depth will also consider- 
ably affect the definition of the core, and although the 
dynamic range of this study is large compared to other 
studies of human body habitats, it is possible that 
microbes we define as 'non-core' still persist at very low 
abundance, perhaps in a manner analogous to a seed 
bank. Determining the role, if any, of these low-abun- 
dance microbes in responding to changes in diet, phy- 
siological status, and so on will be a fascinating 
challenge for future studies. 

Taken together, our observations paint a picture of a 
human microbiome that is highly variable over time 
(Figures 2 and 3; Additional files 2, 3, 4, 5, 6, 7, 8, 9, 
10, 11, 12, 13, 14, 15), but from which stable patterns 
of similarities and differences among body habitats and 
individuals emerge (Figure 1). This temporal variation 
may arise from extrinsic factors, such as exposure to 
different types of foods, medications, or physical envir- 
onments (for example, due to travel), or from intrinsic 
factors, such as the adaptive immune system. Under- 
standing the influence of these factors on an indivi- 
dual's microbiota is an additional challenge for future 
work. 

Although some evidence was found for a high-abun- 
dance core microbiome at the 97% OTU level within 
body sites over time, particularly in the mouth and gut, 
this core appears to represent less than 10% of the total 
OTUs when defining the minimum number of samples 
to be considered part of the core as 90% of the samples. 
This observation does not, however, extend across body 
sites, suggesting that a body-wide core temporal micro- 
biome that bridges the habitat types surveyed here does 
not exist. 



The innovations in sequencing, cloud computing, and 
visualization applied here, together with advances in 
robust and inexpensive microfluidic sample preparation 
techniques, support the development and democratiza- 
tion of inexpensive, informative, and personalized 
microbiome-based phenotyping. Specifically, the ability 
to collect, process and analyze thousands of sequences 
using increasingly available sequencing technologies and 
using commercial computing infrastructure will make 
the ability to trace changes in the microbiome within 
each individual associated with drug administration, dis- 
ease states, and environmental exposures routine. 
Because of the immense subject-to-subject variability in 
the microbiome, studies examining temporal variability, 
which give a view of dynamics beyond the static pictures 
previously available, have the potential to transform our 
understanding of what is 'normal' in the human body, 
and, perhaps, to develop predictive models for the 
effects of clinical interventions. 

Materials and methods 

Sample preparation and sequencing 

Sample collection and DNA isolation were performed as 
described in Costello et al. [6]; and PCR, sequencing, 
and quality filtering of reads were performed as 
described in Caporaso et al. [14]. Samples were not col- 
lected on days 422 through 437. 

To facilitate massively parallel sequencing (1,967 sam- 
ples), barcodes were reused across six lanes in a single 
Illumina GAIIx, with 374, 372, 364, 271, 265, and 323 
samples in lanes 1 through 6, respectively (differing 
from Caporaso et al. [14], where samples were pooled 
and run over seven lanes). Sixteen samples were ulti- 
mately excluded from the analysis as fourteen samples 
were identified as potentially mislabeled (discussed 
below), and the barcodes for two samples were not 
found in the sequencing output, likely indicating a pro- 
blem with amplification for those two samples. 

Data analysis 

To directly compare these M3/F4 time series samples 
with the samples presented in Costello et al. [6], which 
sequenced a different variable region (V2) using a differ- 
ent technology (454 FLX), a reference-based OTU pick- 
ing protocol [18] was applied. After demultiplexing and 
quality filtering sequences, 97% OTUs were picked 
against the Greengenes database [19] (pre-filtered at 
97% identity) using uclust [20]. Reads were assigned to 
OTUs based on their best match to a Greengenes 
sequence, and reads that did not match a Greengenes 
sequence at 97% or greater sequence identity were dis- 
carded. The Greengenes taxonomy associated with the 
best match in Greengenes was assigned to each OTU, 
and the Greengenes tree was used for phylogenetic 
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diversity calculations. These steps and subsequent data 
analysis were performed using Quantitative Insights Into 
Microbial Ecology (QIIME) on AWS. 

Identifying mislabeled samples 

To identify potentially mislabeled samples, we used the 
random forests classifier [21]. A 2,000-tree forest was 
trained on the OTU x Sample Abundance matrix after 
evenly sampling to 500 sequences per sample and 
removing OTUs present in less than 1% of samples. The 
posterior probability that a given sample came from 
each of the body habitats (gut, oral cavity, skin) was esti- 
mated using only those trees in the forest that did not 
contain that sample in their training sets, to avoid over- 
fitting. The classifier considers samples to be mislabeled 
when their alleged environment labels have a low pos- 
terior probability (<60%). Fourteen such samples were 
identified, and these samples were removed from all 
analyses. 

Core microbiome calculation 

The temporal core microbiome across body sites and 
individuals (Figure 2) was computed by varying the 
minimum number of samples in which an OTU must 
be observed to be considered part of the core micro- 
biome, and then determining the number and fraction 
of total OTUs observed in each site (or combination of 
sites) that are part of the core. To facilitate direct com- 
parison across sample types that contained different 
numbers of observations (for example, M3 (all) versus 
M3 gut), we randomly subsampled to exactly 130 obser- 
vations per sample type, corresponding to the sample 
type for which we had the fewest observations. 

Community membership calculations 

The number of consecutive timepoints containing an 
OTU (Figure 3) was calculated as the maximum number 
of consecutive timepoints where an OTU was observed, 
allowing a zero count at a single timepoint to be consid- 
ered part of a continuous stretch of non-zero counts if 
both adjacent timepoints had a non-zero count. This 
controls for sampling error as, for example, a long con- 
tiguous stretch of non-zero counts for an OTU inter- 
rupted by a single zero count for that OTU would likely 
indicate a bad sample, rather than a biologically relevant 
fact about that OTU in relation to the community. Per- 
sistent taxa were defined as those observed in 20% or 
more of the timepoints, but with at least 90% of those 
observations being consecutive (that is, they appear and 
remain present). Transient taxa were defined as those 
observed in at least 60% of the samples, but with at 
most 75% of those observations being consecutive (that 
is, they appear and disappear from the community 
frequently). 



Animated microbial community dynamics 

Animations were created in inVUE [22] based on the 
principal coordinate data presented in Figure la, b. 
inVUE files can be created in QIIME from the principal 
coordinate matrix and associated metadata file. After 
installing and opening inVUE, the user can run, pause, 
and stop the animations associated with different meta- 
data categories. 

Data availability 

All sequence data and sample metadata are publicly 
available under the 'Moving Pictures of the Human 
Microbiome' project [MG-RAST:4457768.3-4459735.3]. 

Additional material 



Additional file 1: Comparison of beta diversity results for 331 
samples sequenced on both 454 and lllumina. Procrustes plot 
comparing principal coordinates of unweighted UniFrac distances. Lines 
connect paired samples sequences on 454 (white tip of line) and 
lllumina (red tip of line). The lllumina samples were evenly sampled to 
5,000 sequences per sample and the 454 samples were evenly sampled 
to 500 sequences per sample. 

Additional file 2: Animation tracing change in position in PCI and 
PC2 with time for all body sites across both individuals. The view 
presented in this video is directly comparable with Figure 1a. 
Background colors correspond to Figure 1a. The M3 time series is shown 
as a red trace (with left palm in orange), and the F4 time series is shown 
as a blue trace (with left palm in white). 

Additional file 3: Animation tracing change in position in PCI and 
PC2 with time for all body sites across M3. The view presented in this 
video is directly comparable with Figure 1a. Background colors 
correspond to Figure 1a. The M3 time series is shown as a red trace 
(with left palm in orange). 

Additional file 4: Animation tracing change in position in PCI and 
PC2 with time for all body sites across F4. The view presented in this 
video is directly comparable with Figure 1a. Background colors 
correspond to Figure 1a. The F4 time series is shown as a blue trace 
(with left palm in white). 

Additional file 5: Animation tracing change in position in PC1, PC2 
and PC3 with time for all body sites across both individuals. 

Background colors correspond to Figure 1a. The M3 time series is shown 
as a red trace (with left palm in orange), and the F4 time series is shown 
as a blue trace (with left palm in white). 

Additional file 6: Animation tracing change in position in PCI, PC2 
and PC3 with time for all body sites across M3. Background colors 
correspond to Figure 1a. The M3 time series is shown as a red trace 
(with left palm in orange). 

Additional file 7: Animation tracing change in position in PC1, PC2 
and PC3 with time for all body sites across F4. Background colors 
correspond to Figure 1a. The F4 time series is shown as a blue trace 
(with left palm in white). 

Additional file 8: Temporal variation in phylum, class, order, family, 
and genus abundances (M3 gut). The x-axis scale differs between M3 
and F4 plots. 

Additional file 9: Temporal variation in phylum, class, order, family, 
and genus abundances (F4 gut). The x-axis scale differs between M3 
and F4 plots. 

Additional file 10: Temporal variation in phylum, class, order, 
family, and genus abundances (M3 tongue). The x-axis scale differs 
between M3 and F4 plots. 
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Additional file 11: Temporal variation in phylum, class, order, 
family, and genus abundances (F4 tongue). The x axis scale differs 
between M3 and F4 plots. 

Additional file 12: Temporal variation in phylum, class, order, 
family, and genus abundances (M3 left palm). The x axis scale differs 
between M3 and F4 plots. 

Additional file 13: Temporal variation in phylum, class, order, 
family, and genus abundances (F4 left palm). The x axis scale differs 
between M3 and F4 plots. 

Additional file 14: Temporal variation in phylum, class, order, 
family, and genus abundances (M3 right palm). The x axis scale 
differs between M3 and F4 plots. 

Additional file 15: Temporal variation in phylum, class, order, 
family, and genus abundances (F4 right palm). The x axis scale differs 
between M3 and F4 plots. 

Additional file 16: Taxonomic summary of the persistent and 
transient OTUs for all individual, body site pairs. 

Additional file 17: Detailed results of persistent versus transient 
community compositions when compared on 454 time series and 
jackknife analysis. 
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